# Import data
df <- readr::read_rds('./data/Juliane_raw_FEST.rds')
# Quick look
dim(df)
## [1] 2735 11
head(df)
## # A tibble: 6 x 11
## # Groups: id [1]
## id block_number trial_no trial_type perc_trial_type CS US
## <dbl> <chr> <int> <chr> <chr> <chr> <chr>
## 1 8.00 baseline 1 CS_only CS_only CSplus <NA>
## 2 8.00 baseline 2 CSminus_only CSminus_only CSminus <NA>
## 3 8.00 baseline 3 CS_only CS_only CSplus <NA>
## 4 8.00 baseline 4 CS_only CS_only CSplus <NA>
## 5 8.00 baseline 5 CSminus_only CSminus_only CSminus <NA>
## 6 8.00 baseline 6 CSminus_only CSminus_only CSminus <NA>
## # ... with 4 more variables: rating <int>, rating_time <int>,
## # button_correct <chr>, button_time <int>
tail(df)
## # A tibble: 6 x 11
## # Groups: id [1]
## id block_number trial_no trial_type perc_trial_type CS US
## <dbl> <chr> <int> <chr> <chr> <chr> <chr>
## 1 30.0 test3 137 CSminus_USlow CSminus_USlow CSminus low
## 2 30.0 test3 138 CSminus_only CSminus_only CSminus <NA>
## 3 30.0 test3 139 CSplus_UShigh CSplus_UShigh CSplus high
## 4 30.0 test3 140 CSminus_only CSminus_only CSminus <NA>
## 5 30.0 test3 141 CSminus_USlow CSminus_USlow CSminus low
## 6 30.0 test3 142 CSplus_UShigh CSplus_UShigh CSplus high
## # ... with 4 more variables: rating <int>, rating_time <int>,
## # button_correct <chr>, button_time <int>
glimpse(df)
## Observations: 2,735
## Variables: 11
## $ id <dbl> 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8...
## $ block_number <chr> "baseline", "baseline", "baseline", "baseline"...
## $ trial_no <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14,...
## $ trial_type <chr> "CS_only", "CSminus_only", "CS_only", "CS_only...
## $ perc_trial_type <chr> "CS_only", "CSminus_only", "CS_only", "CS_only...
## $ CS <chr> "CSplus", "CSminus", "CSplus", "CSplus", "CSmi...
## $ US <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, "low",...
## $ rating <int> 0, -7, -7, -24, -32, -39, -17, -33, -40, -33, ...
## $ rating_time <int> 4824, 4342, 5833, 4173, 2689, 3275, 3460, 3201...
## $ button_correct <chr> "correct", "correct", "correct", "correct", "c...
## $ button_time <int> 8407, 5485, 7576, 5572, 3792, 4610, 4875, 4152...
Ratings on the FEST are plotted for each participant and then for the whole sample. First we plot by trial type delivered, then by trial type perceived.
# Make plots
ggplot(data = df) +
aes(x = trial_no,
y = rating,
colour = trial_type) +
geom_point() +
scale_colour_brewer(palette = "Set1") +
geom_hline(yintercept = 0, linetype = 2) +
geom_vline(xintercept = 10.5, linetype = 3) +
geom_vline(xintercept = 34.5, linetype = 3) +
facet_wrap(~ id) +
labs(title = "Ratings over course of procedure",
subtitle = "Faceted by ID")
df_test <- df %>% filter(trial_type == "CSplus_UStest" | trial_type == "CSminus_UStest")
ggplot(data = df_test) +
aes(x = trial_no,
y = rating,
colour = trial_type) +
geom_point() +
scale_colour_brewer(palette = "Set1") +
geom_hline(yintercept = 0, linetype = 2) +
geom_vline(xintercept = 70.5, linetype = 3) +
geom_vline(xintercept = 34.5, linetype = 3) +
geom_vline(xintercept = 105.5, linetype = 3) +
facet_wrap(~ id) +
labs(title = "Ratings over course of test phases, UStest trials only",
subtitle = "Faceted by ID")
ggplot(data = df_test) +
aes(x = as.factor(trial_type),
y = rating,
colour = trial_type) +
geom_jitter(height = 0) +
geom_boxplot() +
scale_colour_brewer(palette = "Set1") +
facet_wrap(~ id) +
labs(title = "Ratings of test stimuli according to actual trial type",
subtitle = "Faceted by ID")
ggplot(data = df_test) +
aes(x = as.factor(perc_trial_type),
y = rating,
colour = perc_trial_type) +
geom_jitter(height = 0) +
geom_boxplot() +
scale_colour_brewer(palette = "Set1") +
facet_wrap(~ id) +
labs(title = "Ratings of test stimuli according to perceived trial type",
subtitle = "Faceted by ID")
df_test %>% filter(!is.na(trial_type)) %>%
ggplot(data = .) +
aes(x = as.factor(trial_type),
y = rating,
colour = trial_type) +
geom_jitter(height = 0) +
geom_boxplot() +
scale_colour_brewer(palette = "Set1") +
labs(title = "Ratings of test stimuli",
subtitle = "According to actual trial type")
df_test %>% filter(!is.na(perc_trial_type)) %>%
ggplot(data = .) +
aes(x = as.factor(perc_trial_type),
y = rating,
colour = perc_trial_type) +
geom_jitter(height = 0) +
geom_boxplot() +
scale_colour_brewer(palette = "Set1") +
labs(title = "Ratings of test stimuli",
subtitle = "According to perceived trial type")
Here, we plot the prediction interval for the ‘control’ stimulua type (e.g. CSminus/UStest) and depict the trials of the comparator trial type (e.g. CSplus/UStest) that fall within and outside the predictio interval. We do this for the test trials, as delivered and as perceived, and then for the acquisition trials.
# Prediction intervals for test trials as DELIVERED
df_prediction <- df_test %>%
# Filter for CSminus_UStest
filter(trial_type == 'CSminus_UStest' & id != 23) %>%
# Calculate 95% prediction interval
group_by(id) %>%
dplyr::summarise(lower_pi = mean(rating, na.rm = TRUE) - 1.96 * (sd(rating, na.rm = TRUE)),
upper_pi = mean(rating, na.rm = TRUE) + 1.96 * (sd(rating, na.rm = TRUE))) %>%
# Truncate PIs at +50 and -50
mutate(lower_pi = ifelse(lower_pi < -50,
yes = -50,
no = lower_pi),
upper_pi = ifelse(upper_pi > 50,
yes = 50,
no = upper_pi)) %>%
# Join back with df_test
left_join(df_test) %>%
# Is CSplus_UStest within the prediction interval
mutate(CC_effect = ifelse(rating > upper_pi,
yes = 'yes',
no = 'no')) %>%
# Filter for CSplus_UStest
filter(trial_type == 'CSplus_UStest')
ggplot(data = df_prediction) +
aes(x = as.factor(id),
y = rating) +
geom_jitter(aes(colour = CC_effect),
height = 0) +
geom_errorbar(aes(ymin = lower_pi,
ymax = upper_pi)) +
ylim(-50, 50) +
geom_hline(yintercept = 0, linetype = 2) +
scale_colour_brewer(palette = "Set1") +
coord_flip() +
labs(title = 'Prediction interval for CSminus_UStest trial ratings',
subtitle = 'Points are individual CSplus_UStest trial ratings',
y = 'FEST rating',
x = 'Participant ID')
# Prediction intervals for test trials as PERCEIVED
id_19 <- df_test %>% filter(id == 19)
plot(id_19$rating, id_19$trial_no)
df_prediction_perc <- df_test %>%
filter(!is.na(perc_trial_type)) %>%
# Filter for CSminus_UStest
filter(perc_trial_type == 'CSminus_UStest' & id != 23 & id != 19) %>%
# Calculate 95% prediction interval
group_by(id) %>%
summarise(lower_pi = mean(rating, na.rm = TRUE) - 1.96 * (sd(rating, na.rm = TRUE)),
upper_pi = mean(rating, na.rm = TRUE) + 1.96 * (sd(rating, na.rm = TRUE))) %>%
# Truncate PIs at +50 and -50
mutate(lower_pi = ifelse(lower_pi < -50,
yes = -50,
no = lower_pi),
upper_pi = ifelse(upper_pi > 50,
yes = 50,
no = upper_pi)) %>%
# Join back with df_test
left_join(df_test) %>%
# Is CSplus_UStest within the prediction interval
mutate(CC_effect = ifelse(rating > upper_pi,
yes = 'yes',
no = 'no')) %>%
# Filter for CSplus_UStest
filter(perc_trial_type == 'CSplus_UStest')
# Save file df_prediction_perc
write_rds(x = df_prediction_perc,
path = './data/Juliane_raw_as_perc.rds')
ggplot(data = df_prediction_perc) +
aes(x = as.factor(id),
y = rating) +
geom_jitter(aes(colour = CC_effect),
height = 0) +
geom_errorbar(aes(ymin = lower_pi,
ymax = upper_pi)) +
ylim(-50, 50) +
geom_hline(yintercept = 0, linetype = 2) +
scale_colour_brewer(palette = "Set1") +
coord_flip() +
labs(title = 'Prediction interval for ratings of trials perceived as CSminus_UStest ',
subtitle = 'Points are individual CSplus_UStest (percieved) trial ratings',
y = 'FEST rating',
x = 'Participant ID')
# Prediction intervals for all reinforced (CS+/USH or CS-/USL) trials
df_prediction_reinf <- df %>%
# Filter for CSminus_USlow
filter(trial_type == 'CSminus_USlow' & id != 23) %>%
# Calculate 95% prediction interval
group_by(id) %>%
dplyr::summarise(lower_pi = mean(rating, na.rm = TRUE) - 1.96 * (sd(rating, na.rm = TRUE)),
upper_pi = mean(rating, na.rm = TRUE) + 1.96 * (sd(rating, na.rm = TRUE))) %>%
# Truncate PIs at +50 and -50
mutate(lower_pi = ifelse(lower_pi < -50,
yes = -50,
no = lower_pi),
upper_pi = ifelse(upper_pi > 50,
yes = 50,
no = upper_pi)) %>%
# Join back with df_test
left_join(df) %>%
# Is CSplus_UShigh within the prediction interval
mutate(diff_intensities = ifelse(rating > upper_pi,
yes = 'yes',
no = 'no')) %>%
# Filter for CSplus_UShigh
filter(trial_type == 'CSplus_UShigh')
ggplot(data = df_prediction_reinf) +
aes(x = as.factor(id),
y = rating) +
geom_jitter(aes(colour = diff_intensities),
height = 0) +
geom_errorbar(aes(ymin = lower_pi,
ymax = upper_pi)) +
ylim(-50, 50) +
geom_hline(yintercept = 0, linetype = 2) +
scale_colour_brewer(palette = "Set1") +
coord_flip() +
labs(title = 'Prediction interval for CSminus_USlow trial ratings',
subtitle = 'Points are individual CSplus_UShigh trial ratings',
y = 'FEST rating',
x = 'Participant ID')
Calibration was poor. However, participants 11, 15 and 26 have a prediction interval for their CS-/USlow trials that does not cross zero, and they also show no CC effect according to the previous plot.
sessionInfo()
## R version 3.4.3 (2017-11-30)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: OS X El Capitan 10.11.6
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.4/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] bindrcpp_0.2 RColorBrewer_1.1-2 magrittr_1.5
## [4] dplyr_0.7.4 purrr_0.2.4 tidyr_0.7.2
## [7] tibble_1.4.2 ggplot2_2.2.1.9000 tidyverse_1.1.1
## [10] readr_1.1.1
##
## loaded via a namespace (and not attached):
## [1] Rcpp_0.12.15 cellranger_1.1.0 pillar_1.1.0
## [4] compiler_3.4.3 plyr_1.8.4 bindr_0.1
## [7] forcats_0.2.0 tools_3.4.3 digest_0.6.15
## [10] lubridate_1.6.0 jsonlite_1.5 evaluate_0.10.1
## [13] nlme_3.1-131 gtable_0.2.0 lattice_0.20-35
## [16] pkgconfig_2.0.1 rlang_0.1.6.9003 psych_1.7.8
## [19] cli_1.0.0 yaml_2.1.14 parallel_3.4.3
## [22] haven_1.1.0 withr_2.1.1.9000 xml2_1.1.1
## [25] httr_1.3.1 stringr_1.2.0 knitr_1.17
## [28] hms_0.3 rprojroot_1.2 grid_3.4.3
## [31] glue_1.1.1 R6_2.2.2 readxl_1.0.0
## [34] foreign_0.8-69 rmarkdown_1.6 modelr_0.1.1
## [37] reshape2_1.4.3 backports_1.1.1 scales_0.5.0.9000
## [40] htmltools_0.3.6 rvest_0.3.2 assertthat_0.2.0
## [43] mnormt_1.5-5 colorspace_1.3-2 labeling_0.3
## [46] utf8_1.1.3 stringi_1.1.5 lazyeval_0.2.1
## [49] munsell_0.4.3 broom_0.4.2 crayon_1.3.4